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Abstract 


We  deriv  eand  analyze  fast  wideband  algorithms  fdocating  airborne  targets  using  a  passiv  e 
acoustic  sensor  system.  The  system  consists  of  a  distributed  array  of  a  small,  lightweight  acoustic 
vector  sensors,  either  located  in  fixed  positions  on  the  battlefield  or  carried  ty  individual  soldiers. 
These  sensors  measure  the  acoustic  pressure  and  all  three  components  of  particle  velocity  at  a 
single  point  —  the  use  of  such  sensors,  which  have  numerous  ad  van  tages  o  ver  traditional  pressure 
sensors,  has  recently  become  possible  with  the  development  of  a  new  air-based  particle  velocity 
sensor  now  available  commercially  We  deriv  e  a  bearing-only  estimator  based  on  estimation  of 
the  acoustic  intensity  vector  applicable  to  single  vector  sensor.  The  estimator  provides  a  local 
estimate  of  target  bearing  without  the  need  for  inter-sensor  communication,  which  is  useful 
if  the  sensor  is  carried  by  a  battlefield  unit,  for  example..  We  then  develop  a  3-D  location 
estimator  by  combining  the  local  bearing  estimates  at  a  central  processor.  We  also  calculate  an 
optimal  performance  measure  for  the  local  bearing  based  on  the  Cramer-Rao  bound  and  use  it 
to  assess  the  full  potential  of  our  distributed  array  concept  and  the  efficiency  of  our  algorithms. 
Numerical  simulations  are  used  to  show  the  effectiveness  of  our  solutions. 
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1  Introduction 


Acoustic  emissions  from  battlefield  sources  can  provide  an  invaluable  signature  l y  which  to  detect, 
locate,  and  track  hostile  tanks,  trucks,  helicopters,  and  sniper  positions  that  are  camouflaged  or 
have  low  radar  cross  sections.  In  addition,  the  passivity  of  an  acoustic  surveillance  system  allows 
it  to  monitor  the  battlefield  without  giving  aw  ayits  own  presence.  A  nunber  of  researchers  ha  ve 
considered  the  use  of  arra  ysof  pressure  sensors  (ordinary  microphones)  to  perform  this  function 
[1]— [5]  and  have  demonstrated  the  feasibility  of  battlefield  acoustic  localization  and  tracking. 

We  propose  to  consider  the  use  of  ac  oustic  vector  sensors  (AVS’s),  both  individually  and  elements 
of  an  array,  to  perform  the  above  battlefield  surveillance  functions.  These  sensors  measure  the 
(scalar)  acoustic  pressure  and  all  3  components  of  the  acoustic  particle  velocit  yvector  at  a  giv  en 
point.  They  have  been  shown  to  possess  a  number  of  advantages  over  arr^s  ofstandard  pressure 
sensors  that  would  make  them  especially  valuable  on  the  battlefield,  in  both  unattended  and  mobile 
system  contexts.  These  advantages  arise  because  AVS’s  extract  source  information  present  in  the 
structure  of  the  velocit  y  field  that  an  arr§r  of  standard  pressure  sensors  cannot.  F  or  example,  the 
direction  of  the  particle  velocit  yvector  in  a  planewaveor  spherical  w  avocoincides  with  the  line 
from  sensor  to  source.  This  extra  information,  above  and  beyond  a  simple  increase  of  signal-to- 
noise  ratio  (SNR),  is  what  gives  AVS’s  and  AVS  arrays  their  particular  advan  tages  wer  traditional 
pressure-only  systems. 

The  use  of  these  sensors  was  first  considered  from  an  analytical  point  of  view  in  [6],  [7]  for  the 
localization  of  far-held  sources  and  further  studied  in  [8]  and  [9].  A  number  of  estimation  algorithms 
have  also  been  derhed  [10]  and  [11].  F  rom  a  practical  perspective  vector  sensors  have  already  been 
constructed  [12]  and  tested  [14],  [12]— [16] .  Inherent  to  the  current  problem  is  the  presence  of  a 
reflecting  boundary  —  the  ground  —  near  the  sensors.  The  use  of  vector  sensors  near  a  boundary 
has  been  analyzed  in  [17]— [19],  and  they  ha  vebeen  tested  on  a  mock  vessel  hull  [20]  and  at  the 
seabed  [21]. 

Until  now  vector-sensor  research  has  been  directed  to  underwater  situations.  How  ev  er Recently  a 
new  aero-acoustic  probe  called  the  Microhown,  has  become  commercially  available,  from  Microhcwn 
Technologies,  B.V.  [22]  in  the  Netherlands.  This  probe  is  a  micro-machined  sensor  that  measures 
the  differential  resistance  betw  eentw  oheated  cantilevers  in  an  airflow  and  its  output  is  directly 
proportional  to  acoustic  pressure  [23],  [24].  These  sensors  are  currently  available  individually  or 
packaged  with  a  pressure  microphone  to  form  a  single  axis  intensity  probe.  We  propose  combining 
three  orthogonally  oriented  Microflowns  with  a  pressure  sensor  in  a  single  package  to  create  an 
aero-acoustic  vector  sensor.  As  the  probe  is  micro-machined  the  resultant  package  would  be  a  very 
lightw  eigh  sensor  an  amenable  to  being  carried  by  an  individual  soldier. 

In  this  paper  we  develop  decentralized  processing  algorithms  for  locating  airborne  acoustic  sources 
using  a  single  vector  sensor  and  using  a  distributed  arra  y  over  arbitrarily  placed  v  ector  sensors 
located  on  the  surface  of  the  ground.  The  algorithms  are  applicable  to  either  a  fixed  array  of 
unattended  or  a  dynamic  array  of  mobile  sensors  carried  by  individual  soldiers  or  other  battlefield 
units. 

A  decentralized  algorithm  is  one  in  which  some  processing  is  carried  out  locally  in  subarrays  before 
being  combined  centrally  .  In  our  algorithm  the  subarrays  consist  of  individual  vector  sensors  that 
independently  estimate  the  local  bearing  from  their  location  to  the  target.  The  bearingstimate 


is  derived  from  estimates  of  the  three  components  of  the  local  intensity  vector.  An  estimate  is  also 
made  of  its  v  ariabilty.  While  this  information  is  of  use  to  the  soldier  carrying  the  unit  in  its  own 
right  it  can  also  be  transmitted  to  a  central  processor  along  with  its  variability  estimate.  We  then 
propose  to  use  a  weighted  least-squares  (WLS)  techniques  to  combine  the  various  bearing  estimates 
into  a  single  3-D  estimate  of  the  target’s  location.  As  a  variation  w  ealso  propose  a  re-weighted 
LS  method,  to  account  for  the  different  ranges  of  the  sensors  from  the  target.  We  analyze  the 
local  bearing  estimator  b  y  deriving  an  optimal  bound  on  its  mean-square  angular  error  (MS AE) 
for  comparison. 

Any  decentralized  processor  is  suboptimal  since  it  does  not  make  use  of  the  correlations  betw  een 
sensors  at  different  locations,  it  has  a  number  of  practical  advantages.  Firstly  each  sensor  is  small 
enough  to  be  carried  by  an  individual  soldier  which  makes  it  a  very  flexible  system.  There  are  no 
restrictions  on  the  placement  of  sensors  they  may  be  arbitrarily  far  apart,  and  need  not  possess  any 
particular  geometrical  structure.  This  is  an  essential  characteristic  as  the  geometrical  demands  of  an 
array  processing  regime  cannot  be  a  deciding  factor  in  a  unit’s  placement  in  a  hostile  environment. 
F  urthermore,  each  sensor  is  able  to  provide  a  bearing  estimate  of  the  source  to  the  unit  carrying  it 
without  the  need  to  communicate  with  others  sensor  or  the  central  processor  so  not  exposing  itself 
to  detection.  Even  when  communication  is  made  to  the  central  processor,  minimal  data  need  to 
be  sent,  merely  the  bearing  estimate,  an  estimate  of  its  accuracy,  and  the  sensor’s  current  location. 
This  compares  very  favorably  with  a  fully  optimal  design  which  would  require  every  sensor  to 
transmit  every  sample  of  every  component,  thereby  greatly  increasing  the  risk  of  detection  as  well 
as  computational  overhead.  Lastly,  the  algorithms  for  local  bearing  estimation  and  global  position 
estimation  are  both  inherently  wideband  algorithms  and  also  very  computationally  efficient  as  they 
require  no  numerical  optimization.  Thus  estimates  can  be  obtained  very  rapidly,  which  is  almost  a 
de-facto  requirement  of  localizing  airborne  targets,  with  minimal  hardware  costs. 

In  theory  w  ecould  replace  each  vector  sensor  with  a  small  arra  yof  distributed  pressure  sensors 
and  use  the  same  decentralized  processing  methodology,  in  fact  decentralized  algorithms  for  sub¬ 
arrays  of  omni-directional  sensors  have  been  proposed  and  analyzed  for  far-held  [25]  and  near-held 
[26].  How  ever, this  would  require  greater  computational  costs  because  a  numerical  optimization 
procedure  w  ouldbe  required  at  each  array.  In  addition,  the  minimum  frequency  of  interest  in 
airborne  applications  is  50Hz  (see  discussion  in  Section  5)  and  so  eac  h individual  array  w  ouldbe 
incredibly  cumbersome  to  transport  if  it  were  to  be  large  enough  to  attain  any  reasonable  degree  of 
accuracy.  Therefore  the  use  of  vector  sensors  is  essential  to  the  development  of  the  hexible  system 
w  ere  propose. 

In  Section  2  w  epresent  the  mathematical  model  for  the  sensor  measurements,  Section  3  develops 
an  algorithm  to  rapidly  estimate  bearing  using  a  single  vector  sensors.  In  Section  4  w  edevelop 
weighted  and  re-weighted  least-squares  algorithms  for  determining  3-D  aircraft  position  given  the 
bearing  estimates  from  each  sensor,  construct  an  estimator  to  determine  the  weights,  and  giv  e 
an  expression  for  a  lover  bound  on  the  bearing  estimator  from  each  sensor.  Section  5  provides 
numerical  examples  and  Section  6  concludes.  Extensions  are  presented  in  Section  7. 


2  Measurement  Model 


We  assume  that  there  is  a  single  airborne  acoustic  source  radiating  spherical  waves  whose  signal  is 
received  my  m  vector  sensors  at  arbitrary  distinct  locations  on  a  hat  ground  surface. 


If  the  source  is  not  too  close  to  the  boundary,  the  field  due  to  a  point  source  can  be  obtained 
b  y  creating  a  point  image  source,  obtained  b  y  reflecting  the  source  in  the  boundary  and  with 
amplitude  and  phase  determined  by  the  boundary  characteristics,  then  summing  the  fields  from 
the  tw  osources  as  if  the  boundary  w  erenot  present.  Note  that  if  the  source  is  very  close  to  the 
boundary,  ground  waves  and  surface  wi\es  may  exit  [34].  The  resultant  field  may  still  be  obtained 
using  an  image  source  for  locally  reacting  surfaces  but  no  wthe  image  source  must  be  somewhat 
modified  [29].  How  evr,  w  e  shall  not  consider  that  case  in  the  preset  paper. 

Assume,  therefore,  that  the  point  source  is  far  enough  aw  ay  from  ground  that  the  resulting  field  can 
be  regarded  as  arising  from  point  source  radiating  spherically  symmetric  w  a^s  and  simple  point 
image,  therefore  ground  w  a^s  and  surfaces  w  avesdo  not  exist.  This  is  the  case  for  an  airborne 
source.  The  pressure  and  radial  velocity  my  be  described  in  terms  of  the  velocity  potential,  which 
for  a  spherically  symmetric  wave  is 


4>(r,  t)  =  F(t  —  r/c)/r  , 


(2.1) 


where  r  is  the  distance  from  the  source,  c  is  the  sound  speed,  and  F(-)  is  an  arbitrary  differentiable 
function.  The  pressure  and  radial  velocit  y  are  then  [30] 


respectively .  Therefore 
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(2.3) 


(2.4) 

(2.5) 


where  /'(•)  is  the  derivativ  eof  F(-).  Since  the  second  term  in  (2.5)  decreases  as  1/r2  while  the 
first  goes  down  as  1/r,  it  is  negligible  at  distances  of  more  than  a  few  wavelengths.  The  minimum 
frequency  of  interest  in  aircraft  location  is  about  50Hz  corresponding  to  about  6.6m  (see  discussion 
in  Section  5),  so  we  shall  ignore  it. 


Suppose  w  eha  vean  airborne  source  located  at  a  distance  r  and  (unit  length)  bearing  vectors 
from  a  sensor  located  on  the  ground,  which  defines  the  x,y- plane;  the  £  —  axis  is  taken  to  point 
upward.  Construct  an  image  source  by  reflecting  the  source  in  the  x,y-plane  that  also  radiates  the 
spherical  wave  but  with  differed  amplitude  and  phase,  i.e.  if  <f>(r,  t)  represents  the  phasor  notation 
(also  known  as  the  complex  envelope  or  analytic  signal)  of  the  source  field’s  \elocity  potential,  the 
image  source  has  velocity  potential  7 £<f>(r,  £),  where  1Z  is  complex.  Therefore  the  complex  envelopes 
of  the  pressure  and  velocit  y  fields  at  the  sensor  are 

p(t)  =p(t)(l+K)  (2.6) 

v(t)  =  ( 'll  +  7^/iqmag  e)p(t),  (2-7) 

Poc 

where  p(t)  is  the  complex  envelope  of/'(t  —  r/c)/r,  u  is  the  unit  vector  pointing  from  the  sensor  to 
the  (real)  source,  and  U[m3,ge  is  obtained  from  u  by  negating  the  ^-component.  The  presence  of  the 


boundary  imposes  the  condition  —p(t)/vz(t)  =  Z\n  at  £  =  0  [30],  where  vz(t)  is  the  ^-component  of 
v(t),  and  Z\n  is  the  specific  acoustic  impedance  of  the  surface.  Therefore 


_  (1  +  7Z)ppc 

in~  sin^(l-^)’ 

where  ?/>  is  the  elevation  of  the  source  with  respect  to  the  sensor.  Therefore 

^  _  Zjn  -  pqc/  sin  ^ 

Z\yy  +  ppcj  sin?/’’ 


(2.8) 


(2.9) 


The  quantity  1Z  is  known  as  the  reflection  coefficient.  In  general  Z[n  is  a  function,  possibly  quite 
complex,  of  ?/>,  or  equivalently  the  incidence  angle  7  =  pi/ 2  —  ip,  and  frequency.  How  evr,  [34] 
concluded  that  v arious  ground  surfaces  behave  as  if  they  are  locally  reacting.  A  locally  reacting 
surface  is  one  for  which  Z[n  is  independent  of  the  incidence  angle  (but  not  necessarily  frequency)  and 
often  also  arises  in  arc  hitectural  acoustics  with  porous  sound-absorbing  materials  [30].  A  locally 
reacting  surface  may  be  characterized  as  one  in  which  the  sound  disturbance  transmitted  into  the 
low  er  medium  does  not  trarel  along  its  boundary  (actually  this  is  only  strictly  true  for  plane  vaves 
since,  as  mention  above,  ground  wa\es  and  surface  w  a^s  may  exist  when  the  source  is  very  close 
to  the  boundary)  and  therefore  the  normal  velocity  at  each  point  is  completely  determined  by  the 
pressure  at  this  point  [31]. 


It  folio  ws  from  the  abore  development  that  the  measurement  of  a  single  vector  sensor  at  bearing  u 
from  the  source  may  be  written  as  a  complex  four-element  vector 


y(t)  = 


yP(t) 

Vvtt) 


=  hp(t)  +  e(t)  t  =  1,2,... 


(2.10) 


where  yp(t)  is  the  pressure  measurement  yv(t)  contains  the  three  orthogonal  velocity  measurements, 
e(t)  represents  noise,  and  h.  the  steering  vector  is  given  by 


1  +  n 

(1  +  ft)cosc/>cosV’  /911x 

(1  +  1Z)  sin^cos^  ’ 

(1  —  1Z)  sin?/> 

where  </>  is  the  source’s  azimuth  relative  to  the  sensor.  The  expression  (2.11)  for  h  assumes  that 
the  three  velocit  ycomponents  are  aligned  withe  three  coordinate  axes,  or  that  the  orientation 
of  the  sensor  is  known  and  that  the  data  ha  vebeen  rotated  to  achieve  the  same  effect.  It  should 
not  be  difficult  to  design  a  sensor  package  for  which  it  is  easy  talign  the  v  ertical  component.  T  o 
correctly  align  the  horizontal  components  a  compass  would  probably  have  be  included  in  the  sensor 
package,  either  one  that  could  be  aligned  b  y  ey  e  or,  for  sensors  air-dropped  onto  a  battle  field, 
a  digital  compass  whose  output  used  to  rotate  the  data.  Such  an  idea  is  used  in  the  underwater 
v ector  sensors  described  in  [14].  The  measurements  are  again  complex  envelopes.  Note  that  the 
complex  envelope  of  a  bandpass  signal x(t)  is  given  by  x{(t)  +  ixq(t),  where  x{(t)  and  xq(t)  are  the 
in-phase  and  quadrature  component  of  x(t)  respectively  (see  e.g.  [32]). 


We  assume  that  the  signal  and  noise  processes  p(t)  and  e(t)  are  zero- mean  uncorrelated  processes 
with  finite  second  order  moments  and  that 


E  {p(t)p(T)*}  =  a2Jt,T 

E  {e(i)e(r)H}  =  a2I5t,T, 


(2.12) 

(2.13) 


where  <^jT  is  the  Kronecker  delta  function,  I  is  the  identity  matrix,  the  superscript  *  represents 
conjugation  and  superscript  H  complex  conjugation  and  transposition.  In  fact  the  assumption  of 
independent  time  samples  is  not  strictly  necessary  for  the  following  algorithms  to  be  implemented. 
If  there  is  time  correlation,  ho  w  ev  encore  samples  will  be  required  to  achieve  a  given  lev  el  of 
accuracy. 

T  ocalculate  the  performance  measure  for  the  bearing-only  estimator  and  to  implement  the  sim¬ 
ulations  w  ewill  make  the  further  assumption  that  both  processes  are  Gaussian.  If,  for  example, 
the  signal  and  noise  processes  are  band-limited  Gaussian  stochastic  processes,  with  power  spectral 
densities  that  are  symmetrical  about  the  center  frequency,  then  the  complex  envelopes  will  satisfy 
the  above  assumptions  if  they  are  sampled  at  the  Nyquist  rate,  i.e.  twice  the  bandwidth. 

Now  suppose  there  are  m  sensors  located  at  rq, . . .  rm,  and  denote  the  3-D  source  location  vector 
by#.  Letting  p(t)  be  the  pressure  signal  at  the  origin,  the  total  measurements  from  the  array  are 
then 


Vi(t)  =  h(ui)~— - -p(t  -  Ti)  +  ei(t)  i  =  1,...  ,ra,  t  =  1,2,...  (2.14) 

l|0  ”  ri\\ 

where  the  U{  are  the  sensor  to  source  bearing  vectors  and  T{  —  (||#  —  r\\  —  ||0||)/c  is  the  differential 
time  delay  with  respect  to  the  origin.  In  addition  to  the  above  statistical  assumptions  (2.13)  w  e 
also  assume  that  the  noise  processes  are  uncorrelated  from  one  location  to  another  for  the  purposes 
of  simulation,  but  this  is  not  required  for  the  algorithms  to  work.  The  noise  pow  ersat  different 
sensor  locations  may  also  vary  .Note  that  w  edo  not  account  for  relative  Doppler  effects  in  this 
model  or  in  our  simulations,  how  ev  er,  our  estimation  femes  do  not  use  this  information  and  will 
w  ork  just  as  wll  whether  it  is  present  or  not. 

In  the  folio  wing  section  w  ederive  an  algorithm  to  estimate  each  of  the  U{  locally  at  each  vector 
sensor  and  in  Section  4  we  develop  a  methodxtm  bine  these  estimates  to  estimate  0. 


3  Local  Bearing  Estimation 

In  this  section  w  ederive  a  fast  wideband  algorithm  to  estimate  the  bearing  of  the  source  u  from 
a  siigle  v  ector  sensor.  In  the  mobile  array  paradigm,  where  each  sensor  is  carried  by  a  battlefield 
unit,  such  as  a  soldier,  this  estimate  provides  vital  information  the  the  unit  without  the  need  for  any 
communication,  thereby  not  exposing  the  unit  to  detection.  Acoustic  intensity  is  a  vector  quantity 
defined  as  the  product  of  pressure  and  velocity.  Since  the  x  and  y  components  of  the  intensit  y 
vector  are  the  same  for  both  real  and  image  sources  the  acoustic  intensity  vector  is  parallel  to  the 
projection  of  the  source’s  bearing  vector  u  onto  the  x,y-plane.  Therefore  w  ecan  use  an  estimate 
of  the  horizontal  acoustic  intensit  yto  determine  the  azimuth.  Note  that  [7]  used  this  technique 
to  derive  an  estimator  for  the  full  bearing  vector  using  a  vector  sensor  in  free  space,  ho  wver,  it 
cannot  be  used  to  find  the  elevation  when  the  boundary  is  present. 

The  horizontal  component  of  acoustic  intensit  y  measured  ly  a  single  vector  sensor  located  a  distance 
d  from  a  reflecting  boundary  is 

ih(t)  =  yP(t ) 


Vvx  (t) 
VVy(t) ,|  ’ 


(3.15) 


where  the  overbar  indicates  complex  conjugation.  Thus,  for  a  narrowband  or  wideband  source,  as 
long  as  the  noise  at  the  various  sensors  is  mutually  uncorrelated, 


E{Jh(*)} 


a2  |1  +  7£|2  cosip 


COS  (p 

sin  cp 


(3.16) 


Since  this  is  purely  real  w  elet  s  =  N  1  YltLi  Re{/h(^)},  and  by  the  strong  law  of  large  nunbers 
s  E{/h(^)}*  Thus  we  can  estimate  azimuth  from 


.  A 

COS(p 

s  A 

COS  (p 

Uh  = 

si  n(p_ 

—  T7TT7  — ^  Uh  — 

5 

sin  (p 

(3.17) 


Note  that  equation  (3.17)  is  independent  of  1Z ,  and  so  w  e  can  use  this  estiirtsjtcfetermine  the 
azimuth  even  without  knowing  the  local  reflective  properties  of  the  ground.  Since  the  magnitude 
of  the  horizontal  component  of  acoustic  intensity  strongly  depends  on  the  elevation  so  will  the 
accuracy  of  this  azimuthal  estimator.  T  o  be  precise,  ly  adapting  the  analysis  of  £E|etpresen  t 
scenario,  we  can  show  that  its  asymptotic  MSAE  is 

MSAE(«h)  =  lim  N  E(<£  -  4>?  =  0  U  .  •  (3.18) 

v  '  N-> oo  2p2|l  +  7£|4  cos2  ip  v  ' 


T  oobtain  an  estimator  of  the  elevation  is  somewhat  more  complicated.  First  suppose  that  the 
input  impedance  Z[n  is  approximately  constant  across  the  bandwidth  (this  is  the  only  factor  that 
constrains  the  bandwidth  of  our  algorithm  in  an  yw  y)  and  note  that  the  vertical  component  of 

acoustic  intensit  y Iv(t)  =  yp(t)yVz(t)  has  expected  value 

E{Iv(t)}  =  <E2(1  +  ft)(l  -  K)  sin </>  (3.19) 


Using  equation  (3.16)  we  note  that 


a  E{/,(i)}  pc 

||E{/h(*)}||  l  +  11  ZinCOS^’ 


(3.20) 


where  the  second  equality  folio  ws from  (2.9).  The  quantity  on  the  righ  t-handside,  %  say,  is  a 
function  of  pj  alone,  which  we  propose  to  estimate  from  the  statistic 


It  then  remains  to  solve 


££1Re{/1(«)}ir 


X  = 


1  -K{$) 

i  +  n$) 


tanxp 


pc 

Z-m  cos  xp 


(3.21) 


(3.22) 


for  pj.  Thus  we  propose  the  estimate  the  elevation  from 


=  Re 


(3.23) 


We  choose  to  take  absolute  values  before  the  in  verse  cosine,  rather  than  taking  real  values  or 
applying  the  inv  erse  cosine  to  a  complex  number,  because  n  umerical  simulations  sho  w  edhat  it 


resulted  in  slightly  greater  accuracy,  particularly  for  large  N.  The  real  value  of  the  inverse  cosine 
is  tak  en  to  deal  with  (rare)  cases  in  whih  the  argument  is  larger  than  one,  in  these  cases  we  ensure 
that  ty  =  0.  Note  that  the  expression  in  (3.23)  automatically  ensures  that  the  estimated  elevation 
is  not  negative  so  incorporating  the  a-priori  knowledge  than  the  source  lies  above  ground. 

It  can  be  seen  from  the  abov  e  that  the  bearing  is  estimated  via  the  three  components  of  in  ten¬ 
sity.  Therefore,  w  ecould  in  theory  use  three  orthogonally  oriented  single-axis  intensit  yprobes, 
Microflown  Technologies  B.V.  already  packages  its  novel  sensor  with  a  pressure  sensor  to  form  such 
an  intensit  y  probe. In  tensity  probes  are  available  from  other  mamfacturers  such  as  Briiel  and  Kjaer, 
although  those  currently  a  vailable  use^  oclosely  spaced  pressure  sensors,  instead  of  true  velocity 
sensors,  to  determine  the  velocit  y.Such  sensors  are  not  considered  appropriate  for  vector-sensor 
processing  [7]. 


4  Centralized  Location  Estimation 

We  now  consider  the  problem  of  how  to  combine  the  decentralized  estimates  of  the  target  bearings 
to  obtain  an  &4Simate  of  its  location.  Each  sensor  transmits  its  local  estimate  of  the  direction 
from  its  location  to  the  source  iii  say  as  well  as  its  own  location  r*.  The  location  could  be  determined 
from  a  lightw  eigh  GPS  receiver  carried  along  with  the  sensor.  In  practice  both  the  localbearing 
and  position  will  contain  errors,  how  evr,  we  shah  assume  that  the  bearing  estimate  is  the  dominant 
source  of  error  and  ignore  possible  inaccuracies  in  the  location  r\.  Therefore,  a  total  of  five  quantities 
(four  if  ah  sensors  are  at  the  same  altitude)  need  be  transmitted,  three  to  describe  the  location  and 
tw  o  for  the  bearing,  no  matterho  w  long  an  obsermtion  window  is  used.  This  is  a  huge  advantage 
over  a  ceitralized  processing  scheme  in  terms  of  communications  overhead  where  e\ery  single  data 
sample  from  every  single  sensor  must  be  sent  to  a  central  processor. 


4.1  Weighted  Least  Squares 

If  ah  the  ui  w  erewithout  error  the  collection  of  m  lines  passing  through  each  V{  with  direction 
Ui  w  ouldintersect  at  the  true  source  position.  Therefore,  w  ewant  to  choose  the  estimate  of  the 
source’s  location  0  to  be  a  point  that  is  in  sommse  closest  to  ah  these  lines.  We  shah  choose  0 
so  as  to  minimize  a  weighted  sum  of  the  minimum  squared  distances  from  0  to  each  line.  By  doing 
so  we  will  derive  a  closed  form  solution  for  the  location  estimate  so  avoiding  the  need  for  a  complex 
computational  searc  h.  Any  point  along  the  line  defined  b  y  the  iih  sensor’s  location  and  bearing 
estimate  is  defined  by  the  vectors  +  fi{Ui  for  some  /i*.  For  fixed  0  the  point  of  closest  approach 
occurs  when  fii  =  uf  (0  —  r*),  i.e.  fii  is  the  projection  of  the  vector  from  V{  to  0  onto  the  direction 
U{.  Thus  the  weighted  least  squares  (WLS)  estimate  of  0  is 

rri 

0  =  argmin  V"  ||r*  +  uf(0  —  ri)iii  -  0\\2w^  (4.24) 

0  i= i 

where  Wi  is  a  weight  corresponding  to  the  accuracy  of  each  iif.  the  bettethe  accuracy  of  the 
smaller  Wi  should  be.  We  shah  discuss  the  choice  of  weights  fully  in  Section  4.2.  Expanding  (4.24) 
and  rearranging  we  obtain 

rn 

0  =  argminV'  {— 2rJ  (I  —  Uiuj)0  +  6T(I  —  Uiuf)0}  W{ , 

0  i= 1 


(4.25) 


where  we  have  dropped  terms  independeit  of  6.  Note  that  (I  —  iiiiij)  is  the  projection  matrix  onto 
the  plane  orthogonal  to  u\.  Differentiating  with  respect  to  0  and  setting  the  result  equal  to  zero 
gives 

rn 

2^(7  -  Uiuf)(ri  -  0)wi  =  0.  (4.26) 

i— 1 

Hence, 

Aw ,  (4.27) 

where  w  =  [w\, . . .  ,  tcm]T,  W  =  diag{to},  U  =  [u\, . . .  ,  um],  and 

A=[(I  -  «i«i  )n,  •  •  •  ,  (I  ~  UmKn)rm\  •  (4-28) 


e  = 


I  -  uwu 1 


4.2  Choice  of  Weights 

In  general  the  accuracy  of  the  bearing  estimates  ui  will  be  different  from  sensor  to  sensor  due  to  a 
number  of  factors.  There  may  belocal  variations  in  background  noise  level  or  ground  reflectivity, 
signal  strength  will  differ  between  sensors  that  are  different  distances  from  the  source  because  of 
spherical  spreading  loss,  and  the  accuracy  of  the  estimation  algorithm  may  depend  on  the  true 
lx*,  which  differs  between  sensors.  Thus,  it  is  important  for  each  sensor  to  transmit  a  measure  of 
accuracy  only  with  its  bearing  estimate  to  the  central  processor.  A  very  natural  measure  in  this 
situation  is  the  mean-squared  angular  error  (MS AE)  defined  as  E  Sf  where  Si  is  the  angle  betw  een 
Ui  and  its  estimate.  The  analysis  of  our  estimator  is  somewhat  complex  and  no  simple  expression 
is  known  for  the  MSAE.  Although  for  ground  targets,  when  only  azimuth  need  be  estimated,  a 
closed  form  expression  for  the  MSAE  may  be  found. 


Instead  of  the  MSAE  itself,  we  shall  use  a  bound  on  the  MSAE  derived  in  [7],  and  given  by 

MSAEb(w)  =  iV{cos2  V>CRB(</>)  +  CRB(^)},  (4.29) 


where  (f>  and  '</;  are  the  azimuth  and  elevation  of  the  source,  relative  to  the  sensor,  CRB(-)  indicates 
the  Cramer-Rao  bound,  and  N  is  the  nunber  of  snapshots.  It  has  been  shown  that  this  bound 
holds  for  all  N  for  unbiased  estimators  of  a  unit  length  vector  u  and  asymptotically  for  a  muc  h 
larger  class  of  estimators  [27].  In  Appendix  A  we  show  that  the  MSAEb  is  given  by 


MSAEb 


1 

2 P 


cos2  ^  1 

Wdh/d^W2  +  \\dh/dfp\\2  -  \(dhH/dfp)h\2/\\h\\2 


(4.30) 


The  v arious  expression  for  ||fo||2,  \\dh/dip\\2  are  also  given  in  Appendix  A,  from  which  it  w  aybe 
seen  that  MSAEb  is  a  function  of  the  SNR  p  and  the  elevation  ?/>,  but  not  the  azimuth. 


Since  MSAEb  depends  on  the  unknown  quantities  p  and  ^  w  emust  actually  use  an  estimate  of 
MSAEb  by  plugging  in  estimates  of  the  unknowns.  The  bearing  estimator  itself  already  provides 
us  with  and  estimate  of  ^  so  its  remains  to  estimate  p.  If  we  knew  the  maximum-likelihood  (ML) 
estimate  of  u,  say  uMl,  then  w  ecould  use  closed  form  expressions,  see  e.g.  [28],  to  find  the  ML 


estimates  of  a2  and  a2.  Since  we  don’t  know  %l5  we  shall  use  our  actual  estimate  of  ii.  Therefore 
we  have 


tr 


<j?  =  -J—  hH[R-a2I]h, 

llhll4 


(4.31) 

(4.32) 


where  h  =  h(ii )  and  R  is  the  sample  covariance  matrix.  So  our  estimate  of  the  SNR  is  p  —  o2jo2. 
Finally,  we  plug  p  and  the  elevation  estimate  ^  into(4.30)  to  obtain  MSAEb  and  the  weight  that 
w  e  will  use  in  determining#  is  l/MSAEb. 


4.3  Reweighting 

If  a  sensor  is  a  long  w  ayfrom  the  source,  a  small  error  in  its  bearing  estimate  will  have  a  muc  h 
more  dramatic  effect  upon  the  resulting  estimate  of  0  than  if  the  the  sensor  is  very  close  to  the 
source.  Thus,  we  should  also  w  eigh  terms  of  the  squared  error  criterion  more  heavily  for  sensors 
that  are  close  to  the  than  those  far  aw  ay.T  o  be  precise  the  coitiribution  of  the  ith  bearing  estimate 
to  the  squared  error  criterion  is  approximately  where  li  is  the  distance  of  the  source  from  the 
ith  sensor  and  Si  is  the  angular  error  of  U{.  If  w  eknew  the  li  a-priori  w  e would  ideally  like  to 
weight  each  term  by  =  1  /(l2  •  MSAE(u^)). 

Although  we  do  not  know  the  /*,  we  do  have  an  estimate  of  them  after  w  eha  veobtained  the  first 
estimate  of  0  using  the  angular  error  weights  alone,  therefore  w  epropose  a  reweighted  estimator 
constructed  as  follows:  Find  0  from  the  above  WLS  scheme  using  the  angular  error  weights  alone. 
Wi  =  l/MSAE^^).  Using  0,  estimate  the  distances  from  each  sensor  to  the  source  as  U  = 

|| 6  —  r^ll,  then  construct  a  reweighted  estimate  0R  using  the  above  WLS  scheme  but  no  wwith 
Wi  =  l/(ff  •  MSAE^Ui)). 


5  Numerical  Examples 

We  now  illustrate  by  numerical  example  the  performance  of  the  local  bearing  estimator  and  global 
location  estimators  derived  in  the  previous  section.  We  also  suppose  that  the  ground  has  normalized 
input  impedance  Z[n/(poc)  =  11  +  13i,  and  that  it  is  approximately  constant  o\er  this  frequency 
band.  This  corresponds  to  the  value  measured  in  [34]  for  grass-covered  flat  ground  at  225Hz,  the 
mid-point  of  the  band. 


5.1  Local  Bearing  Estimation 

Figure  1  illustrates  the  performance  of  the  local  bearing  estimation  algorithm  via  the  standard  de¬ 
viation  of  its  squared  angular  error,  i.e.  a/MSAE(u)  and  compares  it  with  the  bound  a/MSAE^-m). 
The  scenario  consists  of  a  20dB  source  and  350  snapshots,  and  500  realizations  were  used  for  each 
incidence  angle.  Because  of  the  azimuthal  invarianceof  the  steering  vector  h  of  a  single  vector 
sensor  located  on  the  ground,  all  quantities  are  function  of  elevation  alone. 


Figure  1:  P  erformance, in  terms  of  VMSAE,  of  local  bearing  estimation  algorithm  (solid),  and  bound 
MSAEb  (dashed),  versus  incidence  angle,  for  a  20dB  source  with  350  snapshots.  Also  shown  are  the  averaged 
estimated  MSAEb  (dash-dotted)  plus  three  of  its  empirical  standard  deviations  (dotted).  The  results  were 
obtained  using  500  realizations  at  each  incidence  angle. 


It  can  be  see  that,  at  normal  incidence,  our  algorithms  perfromance  reac  hesthe  MSAEb  but  falls 
off  at  greater  incidence  angles,  underperforming  it  by  about  8°  standard  deviation  by  80°  incidence. 
This  incrseaing  lack  of  optimality  aw  ayfrom  normal  incidence  should  not  have  too  dramatic  an 
effect  in  practice,  how  ev  er^is  it  seems  unlikely  that  the  incidence  angle  would  ever  be  more  than 
about  45°,  even  for  a  low  flying  aircraft  and  a  very  widely  spaced  array.  It  can  be  seen  from  (2.11) 
and  (2.9),  that  \h\  0  as  the  incidence  angle  tends  towards  grazing,  leading  to  very  low  signal 

levels  on  all  component  sensors  at  large  incidence.  This  is  why  all  curves  tend  to  infintity  as  w  e 
appraoch  grazing. 

Not  shown  in  Figure  1  is  the  standard  deviation  of  the  elevation  estimate  as  it  essentially  follows  the 
curve  of  the  MSAE.  Therefore,  we  conclude  that  almost  all  the  angular  error  is  due  to  the  error  in 
estimating  the  elevation  and  that  the  error  in  the  azimuthal  estimates  is  negligible  in  comparison.  T  o 
understand  this  consider  Figure  2.  It  shows  the  signal  gain  resulting  from  reflection  at  the  boundary 
and  the  incoming  signal’s  direction  for  each  sensor  component,  i.e.  the  squared  magnitudes  of  the 
en  tries  of h.  The  pressure  sensor  gain  is  steirtoformly  larger,  because  it  is  omni-directional 
and  is  almost  flat  until  shrply  tapering  off  tow  ardsgrazing  incidence.  The  in-plane  componenet 
gain  (actually  it  is  the  gain  of  the  sum  of  the  in-plane  components,  or  equivalently,  the  gain  of  one 
component  when  u  is  in  the  same  plane  as  it  axis;  this  is  the  quantity  of  interest  because  the  out  put 
of  both  sensors  together  determine  the  azimuth)  is  seen  to  increase  aw  ayfrom  normal,  then  reach 
the  pressure  gain  at  around  80°  before  tapering  off.  It  crosses  the  normal  component  gain  at  about 
5°.  On-the  other  hand,  the  normal  gain,  which  determines  the  ability  to  estimate  elevation,  starts 
off  at  about  —  18dB  and  stays  about  24dB  below  the  pressure  gain  for  all  angles.  The  fact  that  the 
normal  component  gain  is  much  low  er  than  both  pressure  and  in-plane  gain  for  most  angles  explains 
why  the  estimate  of  elevation  contributes  most  of  the  MSAE.  In  addition,  when  the  normal  gain  is 
larger  than  the  in-plane,  very  near  normal  incidence,  large  errors  in  the  azimuthal  estimate  ha  ve 
little  effect  on  the  MSAE  because  of  the  inherent  singularity  in  the  spherical  coordinate  system. 
This  poor  normal  gain  is  mainly  due  to  the  size  of  the  input  impedance,  if  the  input  impedance 
w  esmall,  i.e.  the  surface  were  acoustically  more  pliable,  the  normal  gain  would  improve  relative 


Incidence  Angle  (deg) 

Figure  2:  Signal  gain  resulting  from  reflection  and  sensor  directionality  for  pressure  (solid),  in-plane 
(dashed),  and  normal  (dash-dotted)  components  of  a  vector  sensor  at  a  locally  reacting  boundary.  The 
normalized  input  impedance  is  Z-in  =  11  +  13i. 


to  the  in-plane  and  pressure  gains  resulting  in  better  estimation  of  elevation,  but  poorer  azimuthal 
estimation. 

Also  shown  in  Figure  1  the  average  estimated  bound  MSAEb,  which  will  be  used  as  the  weight 
in  the  global  location  estimator,  plus  3  standard  deviations.  The  estimate  of  the  MSAEb,  is  seen 
to  be  very  accurate  until  about  60°  and  then  begins  to  deviate  slightly  from  tlteue  MSAE  b  and 
starts  to  become  more  variable.  It  is  possible  that  on  any  run  ^  =  0,  with  finite  probablility  if  the 
argument  of  the  inverse  cosine  in  (3.23)  is  greater  than  one.  In  this  case  our  technique  fails  to  yield 
an  estimate  of  MSAEb,  beacuse  it  is  theoretically  infinite  if  really  is  zero  as  \h\  =  0,  so  no  signal 
is  present  on  an  yof  the  sensors.  In  our  simulation,  this  never  occurred  below  66°  incidence,  and 
the  chances  of  it  occuring  rose  to  about  50%  within  a  few  degrees  of  grazing  incidence.  However, 
as  discussed  above,  w  edo  not  expect  that  sue  hlarge  incidences  will  need  to  be  measured  in  this 
application  and  so  this  should  not  be  a  problem  in  practice.  If  it  does  occur  at  one  sensor  in  the 
array,  a  solution  wmld  be  to  use  the  average  of  theMSAEb  obtained  from  the  other  sensors  as  the 
weight  the  failed  sensor. 


5.2  Global  Location  Estimation 

We  no  w  gie  an  example  of  the  performance  of  the  global  location  estimator  for  a  wideband  signal. 
Assume  that  signal  arriving  at  each  sensor  is  a  wideband  stationary  Gaussian  signal  and  that 
it  is  bandpass  filtered  to  have  spectral  support  from  50-400Hz.  Aircraft  noise  above  this  band 
considerably  affected  b  y  propagation  absorption,  while  below  it,  ambient  noise  (predominantly 
caused  by  wind  flo  wingowr  the  sensors)  dominates  [4].  The  signal  is  then  do  wnsampledwith  a 
cen  terfrequency  225Hz  and  the  resulting  complex  envelope  sampled  at  the  Nyquist  frequency  i.e. 
twice  the  bandwidth  700Hz.  We  assume  that  the  signal  and  noise  power  spectra  are  symmetrical 
about  the  center  frequency  so  that  the  resulting  signal  and  noise  samples  are  i.i.d.  complex  Gaussian 
with  variance  equal  to  the  signal  and  noise  pow  ers  respectively. 


Figure  3:  Angular  estimation  performance  in  terms  of  VMSAE  for  the  WLS  and  re- weighted  LS  global 
location  estimators.  The  scenario  is  described  in  the  text. 


Figure  4:  Range  estimation  performance  in  terms  of  VMSRE  for  the  WLS  and  re- weigh  ted  LS  global 
location  estimators.  The  scenario  is  described  in  the  text. 


The  source  is  located  at  (4,5,100.68),  and  there  are  six  sensors  located  on  the  ground  with  x,y 
coordinates  (0,0),  (40,40.52),  (20,18),  (40,-30.53),  (-40,-40.18),  and  (—99.19,0),  all  expressed 
in  meters.  These  locations  were  carefully  chosen,  assuming  a  speed  of  sound  of  330m/s,  to  ensure 
that  the  differential  delays  betw  eenthe  origin  and  each  sensor  are  all  multiples  of  the  sampling 
period  l/700s,  thereby  avoiding  the  need  to  implemeit  fractional  delays  while  generating  data  for 
the  simulation.  The  noise  pow  er  is  assumed  to  be  the  same  at  eah  sensor. 

Figures  3  and  4  show  the  MSAE  and  mean-square  range  error  (MSRE),  defined  by  E(||0||  —  ||0||)2,  of 
the  location  estimated  for  both  the  basic  WLS  estimator  and  the  reweighted  estimator,  versus  the 
SNR  at  the  origin.  A  total  of  350  snapshots  were  used,  corresponding  to  a  half  second  observation 
window  during  which  time  the  source  was  supposed  to  have  been  approximately  stationary.  At 
each  SNR  500  realizations  w  erecomputed.  Unsurprisingly,  both  MSAE  and  MSRE  decrease  as 
SNR  increases,  however,  the  estimated  bearing  is  rather  more  accurate  than  the  estimated  range. 
The  standard  deviation  of  the  angular  error  is  only  3.5°  even  at  OdB  and  falls  below  1°  by  about 
13dB,  whereas  the  standard  deviation  of  the  range  error  is  about  23m  (about  30%  if  its  distance 
from  the  origin)  at  OdB.  It  too  however  is  less  than  5m  by  15dB  and  is  fractions  of  a  meter  at  large 
SNR.  It  is  intuitive  that  a  small  error  in  the  local  bearings  will  cause  a  larger  error  in  the  estimated 
range  than  in  the  estimated  bearing  especia&ythe  source  the  very  far  from  the  sensors.  Hence 
the  somewhat  smaller  angular  errors,  relative  to  range  errors,  seen  in  the  simulation.  Indeed,  it  is 
known  that  the  CRB  on  range  for  a  passive  sensor  array  is  known  to  increase  as  the  fourth  pow  er 
of  the  range  [35]. 

The  reweighted  estimator  is  seen  to  provide  a  slight  decrease  in  both  the  angular  and  range  error, 
the  improvement  is  sligh  tly  greaterfor  the  angle,  at  all  SNR’s.  Given  its  lo  w  computation  cost  it 
is  probably  worth  doing,  even  though  the  gain  is  not  large. 


6  Conclusion 

We  developed  a  fast,  wideband  decentralized  processing  scheme  for  3-D  source  localization  of  air¬ 
borne  targets  using  a  distributed  array  of  aeroacoustic  vector  sensors.  The  procedure  can  be  applied 
to  a  flexible  array,  in  which  each  soldier  or  battlefield  unit  carries  a  lightweight  sensor  package,  or 
to  a  static  array  of  unattended  sensors,  deployed  in  an  air-drop  ever  hostile  territory  perhaps.  The 
algorithm  requires  minimal  communication  betw  een  ea<h  sensor  and  the  base  so  reducing  the  likeli¬ 
hood  of  detection  in  a  hostile  environment.  The  algorithm  proceeds  in  tw  o  stages,  firstly  eah  sensor 
locally  estimates  the  bearing  to  the  source,  using  the  b y calculating  intensities.  These  estimates 
are  then  combined  using  a  WLS  procedure  at  a  central  processor  to  determine  to  3-D  location. 
We  sho  wd  that  performance  could  be  slightly  improved  b  y  using  a  reweighted  procedure  at  this 
second  stage.  We  also  calculated  a  bound  on  the  MSAE  of  the  bearing  estimate  for  each  vector 
sensor  and  sho  wd  that  for  the  most  important  angles  near  normal  incidence  our  algorithm  w  as 
not  too  far  from  optimal. 


7  Extensions 

There  are  a  number  of  direction  in  which  we  will  develop  the  work  presented  in  this  paper,  including 
calculation  of  the  CRB  for  the  full  array  to  pro  vide  a  benchmark  for  comparison  of  the  global 


position  estimator  and  deriving  alternative  decentralized  estimators  based  on  beamforming.  Two 
different  schemes  can  be  considered,  one  in  which  beamforming  is  done  at  each  sensor  to  estimate 
the  local  bearings  and  another  in  which  the  central  processor  is  told  the  covariance  matrix  of  each 
sensor  location,  but  does  not  know  the  cross  covariances  betw  een  locations,  and  maximizes  the  sum 
of  the  individual  spectra.  The  former  requires  more  computational  pow  erto  be  pack  ed  in  to  eac  h 
package  but  each  local  estimate  should  be  a  little  more  accurate  than  the  current  scheme.  The 
latter  requires  minimal  computing  power  at  each  location,  but  requires  transmission  of  more  data, 
namely  the  whole  4x4  co  variancematrix  instead  of  the  bearing.  We  will  also  extend  our  model 
to  account  for  differential  Doppler  shifts  and  develop  tracking  algorithms.  Finally,  we  will  extend 
our  model  to  include  ground  and  surface  w  avdhat  occur  when  the  source  is  close  to  the  ground, 
in  or  to  make  it  applicable  to  ground  targets. 


Appendix  A 


In  this  appendix  we  derive  the  MSAEb  for  a  single  vector  sensor  located  on  the  ground.  The  MSAEb 
defined  by  (4.29)  is  a  function  of  the  CRB  on  azimuth  and  elevation.  Our  measurement  model  for 
a  single  sensor  is  a  member  of  the  class  of  models  whose  CRB  on  the  direction  parameters  is  given 
by  Theorem  3.1  of  [7].  Note  that  this  theorem  only  requires  inversionof  a  2  x  2  matrix  (when 
applied  to  our  case)  even  though  there  are  four  unknowns  as  it  uses  a  tec  hnique  of  concentrating 
the  CRB  with  respect  to  nuisance  parameters  developed  in  [33].  F  or  the  curreit  model  we  obtain 


CRB(0,VO 


(A.l) 


where  D  =  [dh/d(p,  dh/dip]. 


Now,  from  the  definition  of  the  vector-sensor  steering  vector  h  (see  (2.11)),  noting  that  the  reflection 
coefficient  1Z  is  a  function  of  ip  but  not  </>,  we  have  that 
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where  1Z!  is  the  derivative  of  7 Z  with  respect  to  'ip.  Now  (dhH / d(f>)h  =  (dh11  / dcj.))  ( dhdij) )  =  0  follows 
from  (A.2)  and  (A.3),  therefore 
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and  the  result  (4.30)  follows.  F  rom  (2.11),  (A. 2),  and  (A. 3)  ve  also  obtain  the  following  expressions 
for  the  terms  involv  ed  in  MSA® 
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=  \K'{(1  +  1Z)(1  +  cos2  xp)  —  (1  —  IZ)  sin2  p)}  +  (|1  —  1Z\2  —  |1  +  1Z\2)  cos?/’  sin?/?  2  ,  (A. 8) 


where  VJ  indicates  the  complex  conjugate  of  a  VJ . 
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